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INTRODUCTION 


This  is  the  third  in  a  series  of  four  reports  whose  overall  purpose  is  to  describe  an 
adaptive  finite  element  method  (AFEM)  for  solving  systems  of  parabolic  partial  differential 
equations.  In  particular,  AFEM  attempts  to  find  a  ni 
system  of  the  form 

=  [Dix,t,u)u 

subject  to  the  initial  conditions 

afeO)  =  u%x)  , 

and  linear  separated  boundary  conditions 

A\t)u(a,t)  + 

The  variables  x  and  t  represent  spatial  and  temporal  coordinates  and  denote  partial 
differentiation  when  they  are  used  as  subscripts;  u,  f,  u°,  gf,  and  are  M-vectors;  and  D,  A\  A’’, 
and  B  are  MxM  matrices. 

The  problem  is  assumed  to  be  well-posed  and  parabolic;  thus,  e.g.,  D(x,t,u)  is  positive 
definite.  The  rows  of  B  and  B  are  restricted  to  be  either  entirely  zero  or  a  row  of  the  MxM 
identity  matrix.  When  the  row  of  J?'  or  J?"  is  identically  zero,  then  Aj^  or  A[j  cannot  be  zero, 
respectively,  and  the  boundary  condition  is  a  Dirichlet  (essential)  condition.  Otherwise,  the 
boundary  condition  is  a  Neumann  or  Robbins  (natural)  condition.  The  ultimate  goal  of  AFEM 
is  to  determine  an  approximate  solution  to  Eq.  (1)  to  within  a  user  prescribed  error  tolerance. 

The  adaptive  strategies  utilized  by  AFEM  are:  (1)  error  estimation,  coupled  with  (2)  local 
mesh  refinement  (cf.,  e.g.,  Adjerid  and  Flaherty  (ref  1),  Babuska  and  Dorr  (ref  2),  Babuska, 
Zienkiewicz,  Gago,  and  Oliveira  (ref  3),  Bank  and  Weiser  (ref  4),  Berger  and  Oliger  (ref  5), 
Bieterman  and  Babuska  (refs  6,7),  Moore  and  Flaherty  (ref  8),  Shephard  (ref  9),  Strouboulis  and 
Oden  (ref  10),  Zienkiewicz  and  Zhu  (ref  11)),  and  (3)  mesh  movement  (cf.,  e.g.,  Adjerid  and 
Flaherty  (ref  1),  Arney  and  Flaherty  (ref  12),  Bell  and  Shubin  (ref  13),  Davis  and  Flaherty  (ref 
14),  Dorfi  and  Drury  (ref  15),  Dwyer  (ref  16),  Ewing,  Russell,  and  Wheeler  (ref  17),  Hyman  (ref 
18),  Kansa,  Morgan,  and  Morris  (ref  19),  Miller  and  Miller  (refs  20,21),  Petzold  (ref  22),  Rai  and 
Anderson  (ref  23),  Russell  and  Ren  (ref  24),  Saltzman  and  Brackbill  (ref  25),  Smooke  and 
Koszykowski  (ref  26),  Thompson  (ref  27),  Verwer,  Blom,  Furzeland,  and  Zegeling  (ref  28),  and 
White  (ref  29)). 

The  purpose  of  'this  report  is  to  describe  the  mesh  refining  procedures  employed  by 
AFEM.  Detailed  summaries  of  how  AFEM  implements  its  other  adaptive  strategies  are  found  in 
separate  reports  entitled:  Adaptive  Finite  Element  Method  II:  Error  Estimation  (ref  30)  and 
Aadaptive  Finite  Element  Method  IV:  Mesh  Movement  (ref  31).  Furthermore,  the  report. 
Adaptive  Finite  Element  Method  I:  Solution  Algorithm  and  Computational  Examples  (ref  32), 


americal  solution  of  an  Af-dimensional 


,  a<x<b  ,  f>0. 
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describes  how  all  the  adaptive  algorithms  are  implemented  in  unison  and  contains  results 
demonstrating  the  utility  of  AFEM  as  a  computational  tool.  However,  to  fully  understand  how 
and  why  AFEM  performs  mesh  refinement,  a  brief  description  of  its  error  estimating  procedures 
is  presented  here  since  error  estimation  controls  the  mesh  refinement  (cf.,  Coyle  and  Flaherty 
(ref  30)  for  a  more  detailed  account). 

Two  main  error-controlled  adaptive  finite  element  techniques  for  parabolic  partial 
differential  equations  (PDFs)  have  evolved.  One  technique  is  a  method  of  lines  (MOL)  (cf.,  e.g., 
Adjerid  and  Flaherty  (ref  1)  and  Bieterman  and  Babuska  (refs  6,7))  approach  that  uses  a 
Galerkin  method  to  discretize  the  PDEs  in  space  and  one  of  the  many  excellent  ordinary 
differential  equations  (ODEs)  solvers  (cf.,  e.g..  Gear  (ref  33)  and  Petzold  (ref  34))  to  integrate  in 
time.  Spatial  elements  are  added  and  deleted  in  regions  of  high  and  low  error  when  deemed 
necessary,  and  all  elements  are  advanced  in  time  according  to  the  ODE  solver.  The  second 
technique  is  a  local  refinement  method  (LRM)  (cf.,  e.g.,  Berger  and  Oliger  (ref  5)  and  Moore 
and  Flaherty  (ref  8))  that  uses  Galerkin  approximations  and/or  finite  difference  techniques  in 
both  space  and  time.  Fine  grids  of  space-time  elements  are  added  to  coarser  grids,  and  the 
problem  is  recursively  solved  in  regions  of  high  error. 

The  strategy  developed  here  is  an  MOL  technique  with  features  of  an  LRM  which 
attempts  to  take  advantage  of  the  strengths  of  both  methods.  The  basic  idea  is  to  divide  the 
total  error  into  its  spatial  and  temporal  components  (cf.,  Coyle  and  Flaherty  (ref  30)).  Using 
these  various  error  estimates,  one  can  then  decide  on  the  type  of  refinement  that  is  needed: 
spatial,  temporal,  or  both.  MOL  only  refines  in  space,  while  maintaining  the  temporal  error 
tolerance  at  a  fraction  of  the  spatial  error  (cf.,  e.g.,  Adjerid  and  Flaherty  (ref  1)  and  Bieterman 
and  Babuska  (refs  6,7)).  Developers  of  LRM,  on  the  other  hand,  control  the  total  error  by 
simultaneous  refinement  in  space  and  time  (cf.,  Berger  and  Oliger  (ref  5)).  As  a  result,  both 
methods  can  do  more  work  than  necessary. 

Our  algorithm  is  an  MOL  approach  with  the  added  ability  to  "back  track"  and  solve  the 
problem  again  with  a  finer  discretization.  Back  tracking  is  usually  associated  with  an  LRM, 
although  some  MOL’s  codes  now  include  it  (cf.,  Adjerid  and  Flaherty  (ref  1)).  Spatial 
refinement  is  performed  locally  using  the  spatial  error  estimate  to  determine  the  location  and 
degree  of  refinement.  Temporal  refinement  is  global  using  the  temporal  error  estimate  to 
determine  an  improved  time  step.  The  problem  is  solved  again  whenever  refinement  is  necessary 
by  advancing  the  entire  spatial  solution  in  time  with  the  appropriate  time  step. 

This  strategy,  like  MOL  but  unlike  LRM,  avoids  artificial  internal  boundaries  and 
complicated  data  structures  since  there  is  no  need  to  keep  track  of  layers  upon  layers  of  refined 
space-time  grids  as  well  as  the  solutions  on  those  grids  (cf.,  e.g.,  Berger  and  Oligej  (ref  5)  and 
Moore  and  Flaherty  (ref  8)).  However,  like  LRM,  this  strategy  utilizes  the  fact  that  the  ODEs 
arise  from  discretizing  a  PDE,  and  since  this  strategy  has  the  ability  to  back  track,  it  is  more 
difficult  to  deceive  and  less  sensitive  to  coarse  grid  selection  (cf.,  e.g.,  Adjerid  and  Flaherty  (ref 
1),  Berger  and  Oliger  (ref  5),  and  Moore  and  Flaherty  (ref  8)). 

In  order  to  obtain  an  approximate  solution  and  the  various  error  estimates  mentioned 
above,  AFEM  discretizes  Eq.  (la)  in  space  using  Galerkin’s  method  with  piecewise  linear  finite 
elements.  Temporal  discretization  is  performed  by  the  backward  Euler  method.  A  second 
solution  is  calculated  using  trapezoidal  rule  integration  in  time,  and  the  difference  between  the 
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two  solutions  is  used  to  furnish  an  estimate  of  the  temporal  discretization  error.  A  third  solution 
is  obtained  using  quadratic  finite  elements  and  the  trapezoidal  rule  in  time.  This  solution  is 
higher  order  in  space  and  time  than  the  original  piecewise  linear  finite  element-backward  Euler 
solution.  Hence,  it  can  be  used  to  provide  an  estimate  of  the  total  discretization  error  of  the 
piecewise  linear  finite-element  backward  Euler  solution.  Furthermore,  the  difference  between 
the  piecewise  linear  and  quadratic  solutions  calculated  by  the  trapezoidal  rule  furnishes  an 
estimate  of  the  spatial  discretization  error  (cf.,  Moore  and  Flaherty  (ref  8)  or  Coyle  and  Raherty 
(ref  35)). 

AFEM’s  goal  is  to  estimate  the  discretization  error  per  time  step  in  solutions  obtained  by 
using  piecewise  linear  finite  element  approximations  in  space  and  the  backward  Euler  method  in 
time.  It  seems  most  appropriate  to  gage  errors  in  the  norm 

Ml  =  [(«,»«,)  +  ;  (2) 


however,  other  metrics  may  also  be  used.  An  error  estimate  that  is  global  in  space  and  local  in 
time  may  at  first  seem  unusual,  but  it  is  commonly  used  when  spatial  finite  element 
approximations  are  combined  with  temporal  finite  difference  methods  (cf.,  Thom6e  (ref  36)). 

Let  the  piecewise  linear  finite  element  solutions  at  time  f  obtained  by  using  backward 
Euler  and  trapezoidal  rule  temporal  integration  be  denoted  by  If}  and  Uy^,  respectively.  Likewise, 
let  Uy^  denote  the  piecewise  quadratic  finite  element  solution  at  f  found  by  using  the  trapezoidal 
rule  integration  in  time. 

It  is  known  (cf.,  Strang  and  Fix  (ref  37))  that 

-  U"ii')\\i  =  0((AfV)  +  0(N-^)  .  (3) 

Since 

=  0((Af")')  +  CKN-^)  ,  (4) 


one  should  be  able  to  use  the  difference  between  and  U}  to  estimate  the  error  in  U}  ;  thus, 

1“  -  vTi,  s  i{)r*  -  u"tk  *  1“  -  i>;i,  (5) 

£  \K  -  ♦  0((A«")’)  +  . 

The  main  problem  in  using  ||I7%  -  U}\^  as  an  a  posteriori  estimate  of  ||m  -  is  the 
computational  effort  required  to  obtain  Tliis  cost  can  be  reduced  considerably  by  using  the 
superconvergence  property  of  the  finite  element  method  for  one-dimensional  parabolic  systems. 
Nodal  superconvergence  enables  us  to  approximate  Cy^  as 
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(«) 


where  is  a  "hierarehical"  eorrection  term  obtained  by  trapezoidal  rule  integration  and 
pieeewise  quadratie  finite  elements  whieh  vanish  at  the  nodes. 

As  deseribed  above, 


-  ||i4  +  El  - 


nil 


(7) 


furnishes  an  estimate  to  the  error  \\u  -  of  the  baekward  Euler  solution.  Equation  (7) 
suggests  the  inequality 


i-"  ^  wi  -  +  UK 


1/2  111 


(8) 


The  term  -  C^||j  is  the  difference  between  two  piecewise  linear  solutions  computed  with 
temporal  integration  schemes  of  different  orders  and  can  be  regarded  as  a  measure  of  the 
temporal  discretization  error.  In  a  similar  manner,  can  be  regarded  as  a  measure  of  the 
spatial  discretization  error.  Indeed,  when  the  finite  element  system  is  integrated  exactly  in  time, 
Adjerid  and  Flaherty  (ref  1)  proved  that  |l£li  converges  to  the  exact  spatial  discretization  error 
||u  -  t/||i  as  A  -»  00  for  linear  parabolic  problems. 

The  refinement  algorithm  utilizing  the  above  error  estimates  is  described  in  detail  in  the 
next  section.  Following  that,  we  describe  how  error  estimates  are  used  to  determine  the 
necessary  level  of  refinement.  Finally,  a  summary  of  this  report  is  presented. 


MESH  REFINEMENT  ALGORITHM 

The  error  estimates  discussed  in  the  previous  section  are  used  to  control  a  global  mesh 
refinement  procedure  that  maintains  e"  (cf.,  Eq.  (7))  below  a  specified  tolerance  TOL.  Suppose 
that  a  solution  and  error  estimates  and  have  been  calculated  at  time  f  using  a 
mesh  with  N  elements  and  time  step  Af  ^  Further  suppose  that  <  TOL  and  calculate 
solutions  and  error  estimates  at  time  f  =  +  At"  using  a  mesh  with  N  elements  and  time  step 

At"  =  At"'^.  The  refinement  strategy  consists  of  checking  e"  and  proceeding  as  follows: 

1.  If  e"  <  +  llf'yjli  <  TOL,  continue  to  the  next  time  step. 

2.  If  e"  >  TOL,  >  ^2  TOL  and  -  ^2  TOL,  then  determine  a  new 

At"  and  a  new  N,  and  redo  the  integration.  (The  procedure  for  determining  the 
new  At"  and  N  are  illustrated  in  the  next  section.) 

3.  If  e"  >  TOL,  ||t/%  -  f^lli  <  F2  TOL  and  Wv.h  ^  TOL,  then  determine  a  new 
N  and  redo  the  integration. 
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4.  If  e"  >  TOL,  IV?/^  -  ^  V'2  TOL  and  <  V2  TOL,  then  determine  a  new 

Af  and  redo  the  integration. 

Steps  2  through  4  are  repeated  until  step  1  is  satisfied.  If  the  new  AT  and/or  N  are 
properly  determined,  however,  no  additional  repetition  should  be  necessary.  The  choice  of  V2 
TOL  was  determined  as  the  most  efficient  from  various  numerical  experiments. 

The  main  advantage  of  this  refinement  procedure  is  that  separate  estimates  of  the  spatial 
and  temporal  errors  allow  different  strategies  to  be  used  depending  upon  the  dominant 
component  of  the  error.  Thus,  if  the  spatial  component  of  the  error,  as  measured  by  |lEj^||  is 
large,  then  only  spatial  refinement  is  used  to  reduce  the  total  error.  The  opposite  situation  arises 
when  the  spatial  component  of  the  error  is  small  relative  to  the  temporal  component. 

It  is  important  to  note  that  the  error  estimates  used  in  this  refinement  procedure  are,  at 
best,  only  asymptotically  correct.  Thus,  they  will  not  produce  reliable  estimates  on  coarse  meshes 
or  when  errors  are  large.  However,  these  estimates  should  still  indicate  when  refinement  is 
needed  and  how  to  perform  it  (cf.,  Adjerid  and  Flaherty  (ref  1)  and  Moore  and  Raherty  (ref  8)). 

A  local  refinement  strategy,  such  as  those  considered  in  Berger  and  Oliger  (ref  5)  and 
Moore  and  Flaherty  (ref  8),  is  usually  more  efficient  than  the  local/global  strategy  presented 
herein.  The  plan  here  is  to  combine  refinement  with  a  mesh  moving  method  that  equidistributes 
the  change  of  some  local  error  measure  for  a  single  time  step  on  a  mesh  with  a  fixed  number  of 
finite  elements  (cf.,  Coyle  (ref  38)  and  Coyle  and  Flaherty  (ref  31)).  It  may  be  possible  to  use 
this  local/global  refinement  strategy  in  conjunction  with  such  a  mesh  moving  method  since  the 
change  in  the  local  error  measure  will  be  approximately  the  same  on  every  subinterval. 


REFINEMENT  DETERMINATION 

The  choices  for  the  spatial  and  temporal  refinement  levels  are  made  simply  and 
straightforwardly  and  without  undue  computational  effort.  Yet  these  choices  try  to  take 
advantage  of  the  error  estimators  already  computed  and  the  known  asymptotic  nature  of  the  true 
error. 


The  choice  of  a  new  time  step,  At",  is  simple  with  global  temporal  refinement.  Since  the 
local  temporal  error  of  the  Backward  Euler  Method  is  O(Af^),  it  may  be  redueed  by  a  factor  of  r, 
0  <  r  <  1,  by  reducing  A/  by  a  factor  of  r^.  Let  E(At)  denote  the  temporal  error  for  a  time  step, 
At.  If  the  temporal  error  is  such  that  E(At”'^)  >  V2  TOL,  then  one  needs  to  reduce  the  error  by 
some  factor  less  than 


1/2  TOL 

Y  =  _  j 

£(Af"-^)  ’ 

thus,  a  reasonable  choice  for  r  is  %  y.  This  results  in  a  choice  for  At  such  that 


(9) 
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At"  =  r'^At 


'A  A  #n-l 


(0.75)V2TOL 


^l  E(At"-^) 


(10) 


Let  E(x,t)  denote  the  pointwise  spatial  error  function  at  time  t,  then  our  objective  is 


[/: 


E\x^)dx 


VA 


<  Vi  TOL  . 


(11) 


Rewriting  Eq.  (11)  in  terms  of  the  underlying  spatial  discretization  yields 


^  E\x,t)dx 


I'A 


i=l  '■'■  -1 


<  V2  TOL  . 


(12) 


Equation  (12)  implies  a  finer  discretization  is  needed  in  the  interval  (xj_j,Xj)  if 


E\x,t)dx 


> 


VaTOL^ 

N 


TOL^  . 


(13) 


Furthermore,  since  the  true  spatial  error  on  (Xj,j,Xj)  for  piecewise  linear  elements  is  such  that  (cf., 
Wendroff  (ref39)) 


(14) 


where  Nj  is  the  number  of  elements  into  which  (xj,j,Xj)  is  subdivided,  then  to  reduce  the  spatial 
error  by  a  factor  of  r,  one  should  add  r'^  elements  to  (Xj.j,Xj).  This  leads  to  a  choice  of  TV,  such 
that 


Picking  Nj  such  that 

Nj  =  floor  {r-'^)  +  I  (16) 

will  ensure  that  Nj  satisfies  Eq.  (15). 
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Example  1 


Consider  the  following  problem  for  Burgers’  equation: 

Uf  +  UU^  =  €U^  ,0<x<l,^>0, 


u(x,Qi)  =  sin  Ttjc  ,  m(0,0  =  «(1,0  =  0  , 

where  e  =  5x10  ^  The  solution  of  this  problem  is  a  pulse  that  steepens  as  it  travels  to  the  ri^t 
until  it  forms  a  shock  layer  at  x  =  1.  After  a  time  of  0(e''),  the  pulse  dissipates  and  the  solution 
decays  to  zero  (cf.  Figure  1,  at  the  end  of  this  section). 

This  problem  was  solved  for  a  value  of  TOL  equal  to  0.8  but  with  different  initial  values 
for  N  and  A/.  In  Table  1,  a  summary  of  the  results  is  presented  when  N  and  At  are  initially  100 
and  0.1,  respectively.  Table  2  presents  a  summary  for  when  N  is  initially  20  and  A/  is  inidally 
0.02.  Table  3  summarizes  results  when  N  and  At  are  initially  50  and  0.1.  Mesh  trajectories 
corresponding  to  the  data  presented  in  Tables  1,  2,  and  3  are  shown  in  Figures  2,  3,  and  4, 
respectively  (at  the  end  of  this  section).  The  horizontal  lines  in  these  figures  are  used  to 
distinguish  the  regions  where  the  values  of  At  are  different. 

These  tables  present  the  relevant  refinement  data  for  the  significant  time  levels  of  the 
solution  process.  By  a  "significant  time  level,"  we  mean  a  time  level  at  which  a  refinement  has 
occurred  and  one  immediately  preceding  a  refinement.  The  tables  are  divided  into  four  columns 
labelled  "Significant  Times,"  "N"  (with  the  initial  N  in  parentheses),  "At"  (with  the  initial  At  in 
parentheses),  and  "Error  Estimate,"  respectively.  The  significant  time  levels  are  found  in  the  first 
column.  Reading  across  the  tables,  one  finds  the  values  of  N  and  At  necessary  to  complete  the 
solution  step  to  within  the  given  tolerance.  A  dash  is  used  to  indicate  no  change  in  A  or  At  from 
the  preceding  row.  The  last  column  lists  the  value  of  e"  at  the  successful  completion  of  the  time 
step.  The  last  row  of  each  table  contains  results  for  the  final  time  of  0.8.  In  this  manner,  the 
tables  outline  the  refinement  process  giving  some  indication  as  to  why  refinement  occurred,  as 
well  as  demonstrating  the  success  of  the  process. 

These  results  indicate  that  it  is  possible  to  reduce  the  total  error  by  refining  only  in  space 
or  only  in  time,  and  that  the  error  estimates  e",  ||t/%  -  and  can  be  used  to  detect 
when  these  situations  arise. 


(17a) 
(17b, c,d) 
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Table  1  Numerical  Parameters  For  Solving  Example  1  With  TOL  =  0.8  and  N  and  At 

Initially  100  and  0.1,  Respectively 


Significant 

N 

At 

Error 

Times 

(100) 

(0.1) 

Estimate 

0.1000 

- 

- 

0.3221 

0.1566 

- 

0.0566 

0.2218 

0.2131 

- 

- 

0.5436 

0.2549 

165 

0.0417 

0.3455 

0.8000 

- 

- 

0.2728 

Table  2  Numerical  Parameters  For  Solving  Example  1  With  TOL  =  0.8  and  N  and  At 

Initially  20  and  0.02,  Respectively 


Significant 

N 

At 

Error 

Times 

(20) 

(0.02) 

Estimate 

0.24 

- 

- 

0.6741 

0.26 

46 

- 

0.3710 

0.28 

- 

- 

0.4330 

0.30 

66 

- 

0.3072 

0.38 

- 

- 

0.6210 

0.40 

109 

- 

0.2488 

0.80 

- 

- 

0.1595 
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Table  3  Numerical  Parameters  For  Solving  Example  1  With  TOL  =  0.8  and  N  and  St 

Initially  50  and  0.1,  Respectively 


Significant 

Times 

N 

(50) 

At 

(0.1) 

Error 

Estimate 

0.1000 

- 

- 

0.3263 

0.1573 

62 

0.0573 

0.2227 

0.2147 

- 

- 

0.4766 

0.2566 

81 

0.0419 

0.3913 

0.3405 

- 

- 

0.4662 

0.3824 

118 

- 

0.3448 

0.8000 

- 

- 

0.1644 
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SUMMARY 


The  three  error  estimates  in  the  first  section  were  used  to  control  a  global  refinement 
procedure  that  attempts  to  maintain  a  global  measure  of  the  error  per  time  step  below  a 
prescribed  tolerance.  Refinement  can  be  performed  in  space,  time,  or  both  space  and  time 
depending  on  the  dominant  component  of  the  error  estimate. 

The  refinement  algorithm  consisted  of  first  determining  the  three  error  estimates  and 
then  comparing  them  to  the  given  tolerance.  If  the  total  error  estimate  was  small  enough,  no 
refinement  is  performed.  However,  if  the  total  error  estimate  is  large,  its  spatial  and  temporal 
components  are  checked  to  see  if  one  is  more  dominant.  If  one  component  dominates,  then 
refinement  is  performed  only  in  the  corresponding  dimension.  Otherwise,  refinement  is 
performed  in  both  space  and  time. 

The  magnitude  of  the  refinement  is  determined  automatically,  and  a  means  of  doing  this 
is  presented.  This  determination  is  based  on  a  priori  estimates  of  the  numerical  methods 
employed. 

Finally,  Burgers’  equation  was  solved  for  several  initial  values  for  the  temporal  and  spatial 
discretizations  in  order  to  illustrate  the  refinement  algorithm  in  conjunction  with  the  various 
error  estimates. 

The  results  of  Example  1  provide  an  indication  of  the  utility  of  these  estimates  and  of  this 
refinement  procedure.  In  certain  situations,  only  spatial  and  temporal  refinement  were  needed 
to  keep  the  total  error  within  the  prescribed  tolerance,  and  our  error  estimates  could  be  used  to 
determine  when  these  situations  arise. 

This  is  one  of  the  first  attempts  to  simultaneously  address  spatial  and  temporal  errors 
with  different  refinement  strategies.  Most  researchers  (cf.,  e.g.,  Berger  and  Oliger  (ref  5))  have 
used  binary  refinement  in  space  and  time,  but  did  not  attempt  to  determine  the  dominant 
component  of  the  total  discretization  error.  As  noted,  method  of  lines  techniques  (cf.,  e.g., 
Adjerid  and  Flaherty  (ref  1)  and  Bieterman  and  Babuska  (refs  6,7))  typically  assume  that 
temjwral  integration  is  exact  and  refine  based  on  estimates  of  spatial  errors.  There  is  a  great 
potential  for  techniques  that  utilize  different  spatial  and  temporal  refinement  strategies, 
particularly  with  problems  having  singularities. 
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